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The observed single-handedness of biological amino acids and sugars has long been attributed 
to autocatalysis. However, the stability of homochiral states in deterministic autocatalytic systems 
relies on cross inhibition of the two chiral states, an unlikely scenario for early life self-replicators. 
Here, we present a theory for a stochastic individual-level model of autocatalysis due to early life self¬ 
replicators. Without chiral inhibition, the racemic state is the global attractor of the deterministic 
dynamics, but intrinsic multiplicative noise stabilizes the homochiral states, in both well-mixed and 
spatially-extended systems. We conclude that autocatalysis is a viable mechanism for homochirality, 
without imposing additional nonlinearities such as chiral inhibition. 
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One of the very few universal features of biology is 
homochirality: every naturally occurring amino acid is 
left-handed (L-chiral) while every sugar is right-handed 
(D-chiral) [HU]. Although such unexpected broken sym¬ 
metries are well-known in physics, for example in the 
weak interaction, complete biological homochirality still 
defies explanation. In 1953, Charles Frank suggested 
that homochirality could be a consequence of chemical 
autocatalysis |3], frequently presumed to be the mech¬ 
anism associated with the emergence of early life self¬ 
replicators. Frank introduced a model in which the D 
and L enantiomers of a chiral molecule are autocatalyt- 
ically produced from an achiral molecule A in reactions 
A -b D —2 d and A -b L —>■ 2L, and are consumed in a chi¬ 
ral inhibition reaction, D -b L —?> 2A |4]. The state of this 
system can be described by the chiral order parameter w 
defined as w = {d — 1)/(d-\-1), where d and I are the con¬ 
centrations of D and L. The order parameter oj is zero 
at the racemic state, and ±1 at the homochiral states. 
Frank’s model has three deterministic fixed points of the 
dynamics; the racemic state is an unstable fixed point, 
and the two homochiral states are stable fixed points. 
Starting from almost everywhere in the d-L plane, the 
system converges to one of the homochiral fixed points 


(Fig. [1^. 

In the context of biological homochirality, extensions of 
Frank’s idea have essentially taken two directions. On the 
one hand, the discovery of a synthetic chemical system of 
amino alcohols that amplifies an initial excess of one of 
the chiral states [5| has motivated several auto catalysis- 
based models (see |6] and references therein). On the 
other hand, ribozyme-driven catalyst experiments [7], 
have inspired theories based on polymerization and chi¬ 
ral inhibition that minimize [HHin] or do not include at 
all [TT|[I2 autocatalysis. In contrast, a recent experimen¬ 
tal realization of RNA replication using a novel ribozyme 
shows such efficient autocatalytic behavior that chiral in¬ 
hibition does not arise 1131. Further extensions account¬ 


ing for both intrinsic noise [SKIl] and diffusion USHIHj 
build further upon Frank’s work. 

Regardless of the specific model details, all these 
models share the three-fixed-points paradigm of Frank’s 
model, namely that the time evolution of the chiral order 
parameter w is given by a deterministic equation of the 
form [5] 






( 1 ) 


where the function f(t) is model-dependent. However, 
the homochiral states arise from a nonlinearity which is 
not a property of simple autocatalysis, but, for instance 
in the original Frank’s model, is due to chiral inhibition 


(see Fig. Ib). The sole exception to the three-fixed-points 
model in a variation of Frank’s model is the work of Lente 
m, where purely stochastic chiral symmetry breaking 
occurs, although chiral symmetry breaking is only par¬ 
tial, with uj ^ 0 but |a;| < I. 

The purpose of this Letter is to show that effi¬ 
cient early-life self-replicators can exhibit universal 
homochirality, through a stochastic treatment of Frank’s 
model without requiring nonlinearities such as chiral 
inhibition. In our stochastic treatment, the homochi¬ 
ral states arise not as fixed points of deterministic 
dynamics, but instead are states where the effects of 
chemical number fluctuations (i.e. the multiplicative 
noise |20j ) are minimized. The mathematical mechanism 
proposed here is intrinsically different from that 

of the class of models summarized by Eq. ([^. In the 
following, we propose a model which we analytically 
solve for the spatially uniform case and the case of 
two well-mixed patches coupled by diffusion. We then 
show, using numerical simulations, that the results 
persist in a one-dimensional spatially-extended system. 
We conclude that autocatalysis alone can in principle 
account for universal homochirality in biological systems. 
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FIG. 1: (Color online) (a) Phase portrait of Frank’s model: the racemic state is an unstable fixed point (red dot), 
while the homochiral states are stable fixed points (green dots), (b) If chiral inhibition is replaced by linear decay 
reaction, the ratio of D and L molecules stays constant, (c) Adding even the slightest amount of non-autocatalytic 
production of D and L molecules makes the racemic state (green dot) the global attractor of the dynamics. 


Stochastic model for well-mixed system:- Motivated in 
part by the experimental demonstration of autocataly¬ 
sis without chiral inhibition [13] , we propose the reaction 
scheme below, which is equivalent to a modification of 
Lente’s reaction scheme m through the additional pro¬ 
cess representing the recycling of enantiomers: 


A -I- D 2d, 
A D, 

kd 


A-l-L 

A^L. 

kd 


2L, 


( 2 ) 


Compared to Frank’s model, the chiral inhibition is re¬ 
placed by linear decay reactions which model both recy¬ 
cling and non-autocatalytic production. The rate con¬ 
stants are denoted by fc, with the subscript serving to 
identify the particular reaction. The only deterministic 


fixed point of this model is the racemic state (Fig. Icl. 


This model can be interpreted as a model of the evolution 
of early life where primitive chiral self-replicators can be 
produced randomly through non-autocatalytic processes 
at very low rates; the self-replication is modeled by au¬ 
tocatalysis while the decay reaction is a model for the 
death process. 

We now approximate reaction scheme ([^ by means of 
a stochastic differential equation for the time evolution 
of the chiral order parameter, w, which shows that in the 
regime where autocatalysis is the dominant reaction, the 
functional form of the multiplicative intrinsic noise from 
autocatalytic reactions stabilizes the homochiral states. 
We consider a well-mixed system of volume V and total 
number of molecules N. As shown in the Supplemen¬ 
tary Material (SM), for A^ ^ 1, we obtain the following 
equation for w, defined in the Ito sense PO] : 


do; 

dt 


‘2'knkdV /2fcd, 2', /q\ 


Nka \ N 

where r](t) is normalized Gaussian white noise 


The time-dependent distribution of Equation ([^ can 
be computed exactly [24l [23. The stationary distribu¬ 
tion 1^. 


Ps{uj) = Af [l — with a = 


Vkr, 


(4) 


depends on a single parameter, a, where the normaliza¬ 
tion constant N is given by 




(5) 


Equation is compared in Fig. against Gillespie 
simulations [^| of scheme (i). For a = Q!c = 1, w 
is uniformly distributed. For a ^ ac^ where the non- 
autocatalytic production is the dominant production re¬ 
action, Ps{uj) is peaked around the racemic state, a; = 0. 
For a ^ ccc, where autocatalysis is dominant, Ps{uj) is 
sharply peaked around the homochiral states, uj = ±1. 
The simulations were performed for N = 1000, where the 
analytic theory is expected to be accurate; for smaller 
values of N, the theory is qualitatively correct, but very 
small quantitative deviations are observable compared to 
the simulations. For example, for N ~ 100, ac ~ 1.005. 

The deterministic part of Eq. ^ has one fixed point 
at the racemic state, consistently with the phase por¬ 
trait in Fig. Ic The multiplicative noise in Eq. ^ van¬ 
ishes at homochiral states, and admits its maximum at 
the racemic state. For a ac, where autocatalysis is 
dominant, the amplitude of the noise term in Eq. ([^ is 
much larger than the amplitude of the corresponding de¬ 
terministic term. In this regime, the system ends up at 
homochiral states where the noise vanishes. 

To understand this result physically, note that the 
source of the multiplicative noise is the intrinsic stochas- 
ticity of the autocatalytic reactions. While, on average. 
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FIG. 2: (Color online) Comparison between the 
stationary distribution, Eq. Q, (dashed lines) and 
Gillespie simulations of reactions Q (markers), for 
different values of a. Simulation parameters: N = 10^, 

kfi = = kfi = 1. 

the two autocatalytic reactions do not change the vari¬ 
able w, each time one of the reactions takes place, the 
value of uj changes by a very small discrete amount. As a 
result, over time the value of uj drifts away from its initial 
value. Since the amplitude of the noise term is maximum 
at racemic state and zero at homochiral states, this drift 
stops at one of the homochiral states. The absence of 
the noise from autocatalysis at homochiral states can be 
understood by recognizing that at homochiral states, the 
molecules with only one of two chiral states D and L are 
present, hence only the autocatalytic reaction associated 
with that chiral state has a non-zero rate. This reaction 
produces molecules of the same chirality, keeping the sys¬ 
tem at the same homochiral state without affecting the 
value of UJ, and therefore, the variable uj does not expe¬ 
rience a drift away from the homochiral states due the 
autocatalytic reactions. 

Since the stationary distribution of ui in Eq. Q is only 
dependent on a, the decay reaction rate, has no ef¬ 
fect on the steady state distribution of the system. The 
only role of this reaction is to prevent the A molecules 
from being completely consumed, thus providing a well- 
defined non-equilibrium steady state independent of the 
initial conditions. The parameter a is proportional to 
the ratio of the non-autocatalytic production rate, kn, to 
the self-replication rate, ka- In the evolution of early life, 
when self-replication was a primitive function, ka would 
be small and the value of a would therefore be large; but 
as self-replication became more efficient, the value of ka 
would increase and so a would decrease. Therefore, in 
our model, we expect that life started in a racemic state, 
and it transitioned to complete homochirality through 
the mechanism explained above, after self-replication be¬ 
came efficient (i.e. when a ac)- 

It is important to note that all of the previous 


mechanisms suggested for homochirality rely on assump¬ 
tions that cannot be easily confirmed to hold during 
the emergence of life. However, even if all of such 
mechanisms fail during the origin of life, our mechanism 
guarantees the emergence of homochirality, since it only 
relies on self-replication and death, two processes that 
are inseparable from any living system. 

Stochastic model with spatial extension:- We now turn 
to the study of reaction scheme ([^ generalized to the 
spatially-extended case [IT]- We discretize space into 
a collection of M patches of volume V, indexed by i. 
The geometry of the space is defined by (i) — the set 
of patches that are nearest-neighbor to patch i (e.g., for 
a linear chain, (i) = {i — l,i-|-I}). We indicate the 
molecules of species A in patch i by Ai and similarly for 
the other species. Each patch is well-mixed and reac¬ 
tions (|^ occur within, while molecules can diffuse be¬ 
tween neighboring patches with diffusion rate S. In sum¬ 
mary, the following set of reactions defines the spatial 
model: 


k k 

Ai^=A=^Di, Ai^^A^Li, i = 

kd kd 

Ai -h Dj- >■ 2Dj, Ai -h Lj- !■ 2Lj (®) 

T ■ /-N 

T>i — Dj , Li ^— Lj, j € (^). 


We now derive the following set of coupled stochastic 
differential equation for the time evolution of the chiral 
order parameter uji, of each patch i (see SM) 


dwj 

dt 


2knkdV 

Nka 


uji-\- 5 y {ujj — 


jeii) 


Wi) 



{1 - Ujf)TJ,{t) + 



(7) 


where now N represents the average number of molecules 
per patch, t^^’s are independent normalized Gaussian 
white noises, ^i’s are zero mean Gaussian noise with cor- 



EIG. 3: (Color online) Parameter in the 

two-patch system as a function of the diffusion rate <5. 
Gillespie simulations (markers) are compared against 
Eq. © (solid blue line) and Eq. (|T^ (dashed red line). 
Simulation parameters as in Eig. [2] 
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FIG. 4: (Color online) Gillespie simulation of scheme (j^ for a one-dimensional system of M = 100 patches, starting 
from racemic state and ending with all the patches in the same homochiral state uj = —1. Simulation parameters: 

N = 1000, ka = kd = Ij 5 = 10“^, and = 0. 


relator 


k^{i) 


( 8 ) 


+ “ 2) (j) I 6{t - t'), 


and is equal to one if j G (i) and zero otherwise. 

In order to see how the coupling of well-mixed patches 
affects their approach to homochirality, it is instructive 
to consider the simplest case of two adjacent patches 
(M = 2). In the two-patch model, various scenarios 
can happen: the system may not exhibit homochirality 
(oji ^ W 2 ~ 0); each patch can separately reach homochi¬ 
rality (wi = ±1 and UJ 2 = il); the system exbihits global 
homochirality (wi = UJ 2 = ±1). We first analyze the con¬ 
dition for each patch reaching homochirality using per¬ 
turbation theory, in the case of slow diffusion. The sta¬ 
tionary probability density function of the chiral order 
parameter of a single patch, Qsiuj) is defined by 


around a point in space in which the system can be con¬ 
sidered well-mixed. This can be interpreted as the maxi¬ 
mum volume in which diffusion dominates over the other 
terms acting on the variable of interest (in this case w). 
From Eq. 0. this condition is fulfilled for S ~ 2kdOi/N. 
In the vicinity of the transition a is in order of one, there¬ 
fore the condition becomes 5 ^ kd/N. For 5 3> kd/N, the 
whole system can be considered well-mixed, and we can 
find the critical value of a for each patch, starting from 
ttc = 1) from the well-mixed results, and using as volume 
the volume the whole system, i.e., MV. This indicates 
that in a single patch 

^patch^ ford»0. (12) 


A simple formula that interpolates between these ex¬ 
treme limits, asymptotic to 1/M (with M = 2) for large 
6 and to Eq. (111 for small (5, is 


^patch ^ 


5 + 26* 
25 + 25* ’ 



(13) 


^ + 1 

Qs{uj) = J Qs{uj,uj2)duj2 = J Qsiuji,uj)duji, (9) 

where is the joint probability distribution of 

wi and UJ 2 at steady state from Eq. If 5 ^ kd/N or 
smaller, then (see SM) the stationary distribution reads 

Q,(cc)=Z(l-a;2)“+^-\ (10) 

where 2^ is a normalization constant. This result shows 
that the critical a in a single patch, up to the first order 
correction in 5, is given by 

for5«0. (11) 

Ikd 

We can now turn to the case of high diffusion. Recall 
that the patches are defined as the maximum volume 


Figure shows agreement between measured 

from Gillespie simulations of the two-patch system, and 
the Eq. ( |13[ ). At the parameter regime below the ccc curve 
in Fig. individual patches are homochiral. Also, we 
find that the correlation between the homochiral states 
of the two patches increases with diffusion rate 5 and 
become completely correlated when 5 ~ kd/N or more. 
In this regime the system reaches global homochirality. 

This latter result suggests that in the spatially- 
extended model, when autocatalysis is the dominant re¬ 
action (i.e. a is small enough) and when the diffusion rate 
is in the order of kd/N or larger, all patches converge to 
the same homochiral state. Figure]^ shows the dynamics 
of a Gillespie simulation of a one-dimensional chain of 
100 patches, initializes at the racemic state, in the pure 
autocatalytic limit (fc„ —>■ 0). Very quickly, small islands 
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of different homochirality (blue and red) are formed. Is¬ 
lands of opposite chirality competes against each other, 
until the system reaches global homochirality. Note that 
for 6 ~ kd/N we can treat the diffusion process deter¬ 
ministically by ignoring the last term in Eq. 0. In this 
regime, Eq. Q is the same as the equation describing 
one-dimensional voter model, implying that the transi¬ 
tion to homochirality is in the universality class of com¬ 
pact directed percolation [55]. 

In conclusion, a racemic population of self-replicating 
chiral molecules far from equilibrium, even in the ab¬ 
sence of other nonlinearities that have previously been 
invoked, such as chiral inhibition, transitions to com¬ 
plete homochirality when the efficiency of self-replication 
exceeds a certain threshold. This transition occurs due 
to the drift of the chiral order parameter under the in¬ 
fluence of the intrinsic stochasticity of the autocatalytic 
reactions. The functional form of the multiplicative in¬ 
trinsic noise from autocatalysis directs this drift toward 
one of the homochiral states. Unlike some other mecha¬ 
nisms in the literature, this process does not require an 
initial enantiomeric excess. In our model, the homochiral 
states are not deterministic dynamical fixed points, but 
are instead stabilized by intrinsic noise. Moreover, in the 
spatial extension of our model, we have shown that dif¬ 
fusively coupled autocatalytic systems synchronize their 
final homochiral states, allowing a system solely driven 
by autocatalysis to reach global homochirality. We con¬ 
clude that autocatalysis alone is a viable mechanism for 
homochirality, without the necessity of imposing chiral 
inhibition or other nonlinearities. 

T.B. acknowledges valuable discussions with Elbert 
Branscomb. This material is based upon work supported 
by the National Aeronautics and Space Administration 
through the NASA Astrobiology Institute under Cooper¬ 
ative Agreement No. NNA13AA91A issued through the 
Science Mission Directorate. 
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SUPPLEMENTARY MATERIAL 
Well-Mixed System 

We start from reaction scheme (2) of the main paper 


P{x, t), of the system being in the state x at time t pS] , 
We begin by defining the functions F^s as 


— T(^x\x ^S^')P{x 


(S4) 


A + d 2d, A -I- L 2L, 

.kn A kn . ^ 

A D, A L. 

kd kd 


SO that the master equation can be written as: 


(SI) 


dP{x, t) 

m 


- (^^Fm{x+ ^Sm,t) - Fm{x,t'^ ■ (S5) 

m — l ' ' 


Each reaction changes the system from a state x = 
{xi,X 2 ,X 3 ) = {a,d,l), specified by the concentration of 
molecules A, d, and L, to a state of the form x + V~^Sm 
(for some m € 4}), where Sm is the m’th row of 

the stoichiometry matrix 


/-I 

1 


-1 

0 

1 

1 

-1 

0 

V 1 

0 

-1/ 


(S2) 


The probability per unit time of such transition is given 
by transition rates T{x + V~^Sm\x) obtained from law of 


mass action for reactions (SI): 


T(f - 1 - f) = V{kn + kad)a, 


T{x + ^S 2 \x) = V{kn + kal)a, 

T(a:+1 


(S3) 


The set of rates (S3) is used to write the master equation 
for the time evolution of the probability density function, 


This equation defines the stochastic model and can be 
numerically simulated using the Gillespie algorithm [26]. 

In order to initiate an analytical treatment, we begin 
by expanding the right-hand side of the master equation 
(we follow |3D|). We obtain a Gaussian noise approxi¬ 
mation by truncating the expansion at the second-order, 
thus neglecting terms corresponding to higher moments. 
We arrive at the non-linear Fokker-Planck equation: 


dP _ ^ djH^P) 1 A (B^kP) 

dt ^ dxj 2 ^ dxjdxk 

j=i j,k=i ^ 

where the drift vector H with component Hj reads: 


(S6) 




V 


^)l^)^ 


^d(d -t- ^) — (i{2,kn ka{d -\- /)) 
-kdd + a{kn + kad) 


(S7) 


The symmetric diffusion matrix has the form 


^ ~ Y2 


,)|T) {Sm G S, 


V 


kd{d + 1) + a{2kn + ka{d + 1)) —kdd—a{kn + kad) —kdl — a{kn + kal) 

-kdd - a{kn + kad) kdd + a{kn + kad) 0 

(S8) 


where the symbol G indicates the Kronecker product. We 
now decompose the diffusion matrix to B = GG^. Mul¬ 
tiple choices for G exist ESI, and it is easy to check that 
the following 3x2 matrix satisfies the decomposition: 


mean and correlation 


{Vj{t)Vk{t')) =SjkS{t-t'). 


(Sll) 


^ I aJ a (^kad kji) kdd (i(^ka^ -t- k^) -t- kdl 
G = -aJu {kad + kn) + kdd 0 

' 0 -^aikal + krd + kdl 


Note that since, the Fokker-Planck equation (S6) only 


depends on B and not the particular choice of its decom¬ 
position G, the probability density function of x and its 
time evolution do not depend on G either [50] . 


Equation (S6) is equivalent to the following stochastic 


The number of degrees of freedom in Eq. (SIO) can be 


differential equation (defined in the Ito sense) 
da: 


reduced by noting two facts: (i) the reaction scheme (SI) 
conserves the total number of molecules, meaning that 
the total concentration n = a -|- d -I- 1 is conserved; (ii) 
(SIO) simulations show that the concentration r = d + l settles 
to a Gaussian distribution around its fixed point value 
where ryfc’s {k = 1, 2) are Gaussian white noises with zero r*, allowing us to substitute r{t) —)■ r*. We therefore 


dt 


= p[{x) + G{x)fj{t) 
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change variables in Eq. (SIO) (using Ito’s formula) PO] . 
from 


a\ f ( a-\- d + l \ 

d ^ r = d + l , (S12) 

ij \u;J \{d-l)/{d + l) J 


so that the only dynamics occurs in the chiral order pa¬ 
rameter Lo. In the new variables, we find that n = 0 and, 
by taking the positive solution of r = 0, that is 


* _ V (kgri - kd- 2kn)^ + Skgknn + kgU - kd - 2kn 

2k a 


(S13) 

we substitute r —?> r* in the equation for w, and use the 
rule for summing Gaussian variables (i.e. arji -f br ]2 = 
V where a and b are generic functions [5D]) to 
express the stochastic part of the equation using a single 
noise variable. Expressing the result in terms of the total 
number of molecules N = Vn, for ^ 1, we arrive at 
the following stochastic differential equation for chirality 
order parameter w (equation (2) of the main text): 


ai+di+li+a 2 +d 2 +l 2 is conserved; (ii) simulation shows 
that in long time, the variables ri = di + li, r 2 = d 2 + 12 , 
and A = ni — n 2 settle to Gaussian distributions around 
their fixed point values ri = r 2 = r* and A = 0. We do 
the following change of variables 


/ ai \ 


f i^t\ 


( ai + di + li + a2 + d2 + I 2 ^ 

di 


A 


ai + di + li — a2 — d2 — I 2 

h 

_ y 

r\ 


di + li 

02 


r2 


d2 + I 2 

+ 


UJi 


{di — li )/(di -I- li) 

\h J 


V W2 / 


\ id2 — h) / id2 + h) J 


(S18) 


using Ito’s formula. Now the dynamics only occurs only 
in w = (wi,W 2 )- For large average number of molecules 
per patch TV ^ 1, the resulting Fokker-Planck equation 
for time evolution of the joint probability density function 
of wi and UJ 2 , Q(uj,t), reads 


dQ ^ d{{Ld;)^Q) I " {U,,Q) 

dt ^ duji 2 ^ dujidujj 


dw 

dt 


2knkdV 

Nka 


w -I- 


N 


(I-a;2)7y(t), (S14) 


where rj{t) is Gaussian white noise with zero mean and 
unit variance. The corresponding Fokker-Planck equa¬ 


tion of Eq. (SI4) is an exactly solvable partial differential 
equation with time dependent solution given in [25j . The 
steady state probability distribution of w is given by 




Ps(w) = Af (l — with a = 


with the normalization constant 




(S15) 


(S16) 


Note that the above sums are now over the patches, and 


not over species as in Eq. (S6). The Jacobian matrix L 


L = - 


2kdknV 

Nk„, 


1 0 
0 1 


-1 1 
1 -1 


and the diffusion matrix U reads 
2fcd ( 1 — wi 0 


U = 


N 


1 — ojI 
0 


1-ujI 


N 


2(1 — ijJ\UJ2') + (^2 — ^ 

4“ ^2 — ^ — ^ 1 ^ 2 ')- 


(S20) 


(S21) 


Note that the stochastic differential equation correspond¬ 


ing to Eq. (S19) is the M = 2 case of equation (7) of the 


mam paper. 


Two-patch model 

Starting from reaction scheme (6) of the main paper 
for the spatial extension of our model 

k /c 

Di, Ai Lj, i = 

kd kd 

A, + Di-^2D„ A, + Li^2hi (S17) 

- Dj , - Lj , J G { 1 } , 

for M = 2, we can follow the procedure explained in the 
previous section to obtain a Fokker-Planck equation for 
time evolution of the probability density of system being 
at a state with concentrations ai, di, li, 02 , d 2 , and l 2 - 
Again we can reduce the number of variables using the 
following facts (i) the total concentration nt = ni+n 2 = 


Perturbation theory for the small diffusion case 


Now, we can use perturbation theory for small 6, to 
find the stationary probability density of ui for each patch 
defined as 


P+1 P+1 

Qs{uj) = / Qsiuj,uj2)duj2 = / Qs{uJi,uj)duji, 

(S22) 

where Qs{uji,uj 2 ) is the stationary solution of Eq. ( S19| ). 
For S ^ kd/N or smaller, we can treat the diffusion de¬ 
terministically by ignoring the last term in Eq. (S211. To 
solve for Qg{u)), we begin by rewriting Eq. (S19) as a 
continuity equation, 


atQ-t-V-J = 0, 


(S23) 































which defines the probability current J as 


The last integral can be evaluated using Bayes’ theorem 


J = Lw Q - -V • (UQ). 


(S24) 


By the conservation of probability, at stationary condi¬ 
tions, the total probability flux through each vertical 
section of a;i-a ;2 plane must be zero. That is 


J Js,ldL02 — J ^{'Llo)iQs — -dujj^{UiiQs)^ dL02 


(S25) 


f+i 

+ W2Qs(wi, W2)dW2 = 0. 


/ +J- P + i- 

U}2Qs{‘^lT^2)du!2 = 6 / UJ2Qs(,^2\^l)Qs{^l)du!2 

= SQsiuJi){uj2)u.,=0{6^), 

(S26) 

which is of order S'^ for small 6, since, {oJ 2 )u-i (the ex¬ 
pected value of L 02 given uji) vanishes at zero 5, and 


therefore, of order <5 for small <5. In this regime, Eq. (S25) 


provide us with a differential equation for Qs(w) with the 
solution (equation (10) of the main paper) 


--1 


Qs[uj)=Z{l-u:y 
where the normalization constant Z is given by 


Z = 


r (« + C + 5 ) 


v^r(a+ 2 ^) 


(S27) 


(S28) 




